package com.github.zyf.ioaj.algorithm.swarmintelligence;

import com.github.zyf.ioaj.constant.FileType;
import com.github.zyf.ioaj.export.Report;
import com.github.zyf.ioaj.export.Tracker;
import com.github.zyf.ioaj.algorithm.Input;
import com.github.zyf.ioaj.algorithm.Output;
import com.github.zyf.ioaj.constant.ParamType;
import com.github.zyf.ioaj.annotation.internal.Algorithm;
import com.github.zyf.ioaj.annotation.internal.BasicParam;
import com.github.zyf.ioaj.annotation.internal.TempParam;
import com.github.zyf.ioaj.utils.basedata.ArrayUtils;

/**
 * <h1>Particle Swarm Optimization</h1>
 * <p>粒子群优化算法</p>
 * <p>默认求解最小值，如果需要求解最大值，将目标函数乘以-1即可，最后得到的值再次乘以-1就是得到的最优值（自变量取值不用乘以-1）</p>
 * <p>需要计算则调用{@link #go()}</p>
 * <p>需要跟踪每次迭代则调用{@link #tracker(String, String, String, FileType)}</p>
 * <p>需要生成报告则调用{@link #report(String, String, String, FileType)}</p>
 * <p>由于算法默认采取了自适应速度惯性权重优化以及压缩因子优化，因此速度惯性权重以及学习因子不可自定义，如有需要只能自定义粒子个数与迭代次数</p>
 * <p>粒子群算法基本思想如下：</p>
 * <p>
 * 现在有一群鸟（一只鸟就好比一个粒子）需要觅食，我们对这群鸟作出几个假设：
 * <ol>
 *     <li>所有鸟都不知道食物最丰盛的地方在哪</li>
 *     <li>所有鸟同时共享当前已知最为丰盛的食物所在地方（谁找到哪边最丰盛，那大伙就都知道了）</li>
 *     <li>
 *         每只鸟有一定的“判断能力”。对于其中的某一只鸟它的“选择”（并非主观选择，由随机数控制）有三种：
 *         <ul>
 *             <li>惯性：不做出任何调整，原来往什么方向飞，现在还是什么方向</li>
 *             <li>从众：做出一个便捷从众的调整，向着目前最优的位置飞</li>
 *             <li>自主：做出一个超级自信的调整，向着自己去过的最优位置飞</li>
 *         </ul>
 *     </li>
 * </ol>
 * <p>
 * 具体计算方法如下：
 * <ul>
 *     <li>速度更新：v = w * v + r1 * c1 * (pbest - p) + r2 * c2 * (gbest - p)</li>
 *     <li>位置更新：p = p + v</li>
 *     <li>
 *         速度惯性权重自适应变化：
 *         <ul>
 *             <li>f(p) \ge f_average: w = w_max</li>
 *             <li>f(p) \le f_average: w = w_min + (w_max - w_min) * ((f(p) - f_min) / (f_average - f_min))</li>
 *         </ul>
 *     </li>
 *     <li>压缩因子：c1 = c2 = 2.05, phi = 0.7298</li>
 * </ul>
 * 其中
 * <ul>
 *     <li>v：{@link #v}</li>
 *     <li>p：{@link #p}</li>
 *     <li>w：{@link #w}</li>
 *     <li>c1, c2：{@link #c1}, {@link #c2}</li>
 *     <li>r1, r2：两个随机数</li>
 *     <li>pbest：{@link #pbest}</li>
 *     <li>gbest：{@link #x}</li>
 *     <li>f：{@link #fit}, {@link Input#getFunction()}</li>
 *     <li>f_average：当前适应度均值</li>
 *     <li>f_min：当前适应度最小值</li>
 *     <li>w_min, w_max：预先给定的最小惯性系数和最大惯性系数，一般取0.4和0.9</li>
 * </ul>
 * <p>
 * 在这里，我们以求最小值为例（如果是最大值就反过来）：
 * 适应度越小，说明距离最优解越近，此时更需要局部搜索；
 * 适应度越大，说明距离最优解越远，此时更需要全局搜索。
 * </p>
 *
 * @author zyf
 */
@Algorithm
public class PSO extends SwarmIntelligenceAlgorithm {

    /**
     * 每一个粒子的当前个体最优值，double[n][input().varNum()]
     */
    @TempParam(type = ParamType.DOUBLE_ARRAY_2, desc = "每一个粒子的当前个体最优值")
    private double[][] pbest;

    /**
     * 当前群体最优值，double[input().varNum()]
     */
    @TempParam(type = ParamType.DOUBLE_ARRAY_1, desc = "当前群体最优值")
    private double[] x;

    /**
     * 每一个粒子当前所在位置，double[n][input().varNum()]
     */
    @TempParam(type = ParamType.DOUBLE_ARRAY_2, desc = "每一个粒子当前所在位置")
    private double[][] p;
    /**
     * 每一个粒子当前速度，double[n][input().varNum()]
     */
    @TempParam(type = ParamType.DOUBLE_ARRAY_2, desc = "每一个粒子当前速度")
    private double[][] v;
    /**
     * 每一个粒子当前适应度，double[n]
     */
    @TempParam(type = ParamType.DOUBLE_ARRAY_1, desc = "每一个粒子当前适应度")
    private double[] fit;


    /**
     * 速度惯性权重
     */
    @TempParam(type = ParamType.DOUBLE, desc = "速度惯性权重")
    private double w;
    /**
     * 个体学习因子
     */
    @BasicParam(type = ParamType.DOUBLE, desc = "个体学习因子")
    private double c1;
    /**
     * 群体学习因子
     */
    @BasicParam(type = ParamType.DOUBLE, desc = "群体学习因子")
    private double c2;

    /* 初始化部分 */

    /**
     * <p>默认进行算法优化</p>
     *
     * @param input 输入量 {@link Input}
     */
    public PSO(Input input) {
        super(input);
        init();
    }

    /**
     * 初始化
     */
    private void init() {
        /* 设置基本的参数 */
        // 惯性权重初值（动态变化）
        this.w = 0.9;
        // 个体学习因子（不可更改）
        this.c1 = 2.05;
        // 群体学习因子（不可更改）
        this.c2 = 2.05;

        /* 为每个数组分配空间 */
        // 个体最优解
        this.pbest = new double[getN()][getInput().getVarNum()];
        // 全局最优解
        this.x = new double[getInput().getVarNum()];
        // 粒子当前位置
        this.p = new double[getN()][getInput().getVarNum()];
        // 粒子速度
        this.v = new double[getN()][getInput().getVarNum()];
        // 粒子适应度
        this.fit = new double[getN()];

        int i, j;

        /* 初始化每个粒子的位置，以及个体最优解和全局最优解 */
        for (i = 0; i < getN(); i++) {
            for (j = 0; j < getInput().getVarNum(); j++) {
                p[i][j] = getInput().randomValueInRange(j);
                pbest[i][j] = 0.0;
            }
        }
        x = p[0].clone();

        /* 初始化每个粒子速度 */
        for (i = 0; i < getN(); i++) {
            for (j = 0; j < getInput().getVarNum(); j++) {
                // 粒子的最大速度通常为可行域的10% - 20%
                double max_v = getInput().getRange(j) * 0.2;
                // 控制所有速度在-max_v到max_v之间
                v[i][j] = -max_v + 2 * max_v * Math.random();
            }
        }

        /* 根据初始化值经计算适应度、个体最优解以及全局最优解 */
        for (i = 0; i < getN(); i++) {
            fit[i] = getInput().getFunction().apply(p[i]);
            // 对gbest赋初值
            if (fit[i] < getInput().getFunction().apply(x)) {
                x = p[i].clone();
            }
            // 以第一组计算值作为每个粒子目前的最优值
            pbest[i] = p[i].clone();
        }
    }

    /* 计算部分 */

    /**
     * 每次迭代时扫描一遍
     */
    public void scan() {
        int i, j;
        // 计算每个粒子的平均适应度和最小适应度
        double f_average = ArrayUtils.mean(fit);
        double f_min = ArrayUtils.min(fit);
        // 第二层循环，扫描所有粒子
        for (i = 0; i < getN(); i++) {
            // 第三层循环，扫描所有变量
            for (j = 0; j < getInput().getVarNum(); j++) {
                // 计算自适应权值
                // 大于平均值需要全局搜索
                if (fit[i] > f_average) w = 0.9;
                    // 小于平均值需要局部搜索
                else w = 0.4 + 0.5 * (fit[i] - f_min) / (f_average - f_min);

                // 更新速度
                v[i][j] = w * v[i][j] + c1 * Math.random() * (pbest[i][j] - p[i][j]) + c2 * Math.random() * (x[j] - p[i][j]);
                // 压缩因子
                v[i][j] *= 0.7298;

                // 限制速度范围，避免容易越界
                double max_v = (getInput().getUb()[j] - getInput().getLb()[j]) * 0.2;
                if (v[i][j] < -max_v) v[i][j] = -max_v;
                if (v[i][j] > max_v) v[i][j] = max_v;

                // 更新粒子位置，时间单位为1，故这里省去”*1“
                p[i][j] += v[i][j];

                // 调整位置，避免越界
                if (p[i][j] < getInput().getLb()[j]) p[i][j] = getInput().getLb()[j];
                if (p[i][j] > getInput().getUb()[j]) p[i][j] = getInput().getUb()[j];
            }

            // 计算目标函数值（适应度）
            fit[i] = getInput().getFunction().apply(p[i]);

            // 更新全局最优位置
            if (fit[i] < getInput().getFunction().apply(x)) x = p[i].clone();
            // 更新每个粒子的最优位置
            if (fit[i] < getInput().getFunction().apply(pbest[i])) pbest[i] = p[i].clone();
        }
    }

    /**
     * 生成最优解
     */
    private void generateOptimization() {
        // 获取最优值
        double optimalOutput = getInput().getFunction().apply(x);
        setOutput(new Output(optimalOutput, x));
    }

    /**
     * 运行并计算结果
     */
    public void go() {
        int t;
        for (t = 0; t < getIter(); t++) {
            scan();
        }
        generateOptimization();
    }

    /**
     * 运行并计算结果（按迭代逐步跟踪），如果目的是计算出结果，请使用{@link #go()}
     *
     * @param author   作者
     * @param message  备注
     * @param title    标题
     * @param fileType 文件类型 {@link FileType}
     * @return {@link Tracker}
     */
    @Override
    public Tracker tracker(String author, String message, String title, FileType fileType) {
        int t;
        Tracker tracker = new Tracker(this, author, message, title, fileType);
        tracker.startTracker();
        for (t = 0; t < getIter(); t++) {
            scan();
            tracker.trackOne(t + 1);
        }
        generateOptimization();
        tracker.endTracker();
        return tracker;
    }

    /**
     * 生成报告
     *
     * @param author   作者
     * @param message  备注
     * @param title    标题
     * @param fileType 文件类型 {@link FileType}
     * @return {@link Report}
     */
    @Override
    public Report report(String author, String message, String title, FileType fileType) {
        return new Report(this, author, message, title, fileType);
    }

    /* GETTER & SETTER 部分 */

    public double getW() {
        return w;
    }

    public double getC1() {
        return c1;
    }

    public double getC2() {
        return c2;
    }

    /* GETTER 部分 */

    public double[][] getPbest() {
        return pbest;
    }

    public double[] getX() {
        return x;
    }

    public double[][] getP() {
        return p;
    }

    public double[][] getV() {
        return v;
    }

    public double[] getFit() {
        return fit;
    }
}
